DOCUMENT RESUME 



ED 428 100 



TM 029 501 



AUTHOR 

TITLE 

INSTITUTION 

PUB DATE 
NOTE 

AVAILABLE FROM 

PUB TYPE 
EDRS PRICE 
DESCRIPTORS 



IDENTIFIERS 



Beguin, Anton A. ; Glas, Cees A. W. 

MCMC Estimation of Multidimensional IRT Models. Research 
Report 98*14. 

Twente Univ. , Enschede (Netherlands) . Faculty of Educational 
Science and Technology. 

1998-00*00 
3 Op . 

Faculty of Educational Science and Technology, University of 
Twente, P.O. Box 217, 7500 AE Enschede, The Netherlands. 
Reports - Evaluative (142) 

MF01/PC02 Plus Postage. 

Ability; *Bayesian Statistics; *Estimation (Mathematics) ; 
*Item Response Theory; *Markov Processes; *Monte Carlo 
Methods; Simulation 

Gibbs Sampling; *Multidimensionality (Tests) ; Two Parameter 
Model 



ABSTRACT 



A Bayesian procedure to estimate the three -parameter normal 
ogive model and a generalization to a model with multidimensional ability 
parameters are discussed. The procedure is a generalization of a procedure by 
J. Albert (1992) for estimating the two-parameter normal ogive model. The 
procedure will support multiple samples from multiple populations and 
restrictions on the factor matrix for testing specific hypotheses about the 
ability structure. The technique is illustrated using simulated and real 
data. A Markov chain Monte Carlo (MCMC) procedure is used to sample the 
posterior distributions of interest and needed chains are constructed using 
the Gibbs sampler (A. Gelfand and A. Smith, 1990) . (Contains 8 tables and 47 
references . ) (Author/SLD) 



******************************************************************************** 

* Reproductions supplied by EDRS are the best that can be made * 

* from the original document. * 

******************************************************************************** 



TM029501 



TK 



O 

O 



00 

(N 



s 



MCMC Estimation of 
Multidimensional IRT Models 



Research 

Report 

98-14 



□ 



PERMISSION TO REPRODUCE AND 
DISSEMINATE THIS MATERIAL HAS 
BEEN GRANTED BY 

o. ^-ej 1 S3 



Anton A. Beguin 
Cees A.W. Glas 



TO THE EDUCATIONAL RESOURCES 
INFORMATION CENTER (ERIC) 



originating it P n or organization 



Oocumenu^Tnot ns'cess* s . ,a,e<< in this 
official OERI position 0r p 0 r ;^ y represen * 





faculty of 

EDUCATIONAL SCIENCE 
AND TECHNOLOGY 




Department ofi 1 

Educational Measurement and Data Analysis 



2 




University of Twente 



MCMC estimation of multidimensional IRT models 



Anton A. Beguin 



Cees A.W. Glas 

Department of Educational Measurement and Data Analysis 
Faculty of Educational Science and Technology 
University of Twente, Enschede, the Netherlands 
E-mail A.A.Beguin@edte.utwente.nl 




3 



Abstract 



In this paper, a Bayesian procedure to estimate the 3 parameter normal ogive model 
and a generalization to a model with multidimensional ability parameters will be discussed. The 
procedure is a generalization of a procedure by Albert (1992) for estimating the 2 parameter 
normal ogive model. The procedure will support multiple samples from multiple populations 
and restrictions on the factor matrix for testing specific hypotheses about the ability structure. 
The technique is illustrated using simulated and real data. 

Key words: Bayes estimates. Full Information Factor Analysis, Gibbs sampler, item 
response theory, Markov chain Monte Carlo, multidimensional item response theory, normal 
ogive model. 

Introduction 

Item response theory (IRT) models are stochastic models for the responses of persons 
to items, where the influence of items and persons on the responses are modeled by disjunct 
sets of parameters. In the framework of educational and psychological measurement, the person 
parameter can usually be labeled ability, and in many instances it suffices to assume that ability 
is unidimensional. However, in other instances it may be a priori clear that multiple skills or 
abilities are involved in producing the manifest responses, or the dimensionality of the ability 
structure might not be clear at all. In these situations, multidimensional IRT models can serve 
confirmatory and explorative purposes. The focus of this paper will be dichotomously scored 
items. Multidimensional IRT models for dichotomously scored items were first presented by 
Lord and Novick (1968) and McDonald (1967). These authors use a normal ogive to describe 
the probability of a correct response. McDonald (1967,1997) developed an estimation procedure 
based on an expression for the association between pairs of items and derived from a polynomial 
expansion of the normal ogive. The procedure is implemented in NOHARM (Normal-Ogive 
Harmonic Analysis Robust Method, Fraser, 1988). An alternative using all information in 
the data, and therefore labeled ’’Full Information Factor Analysis” was developed by Bock, 
Gibbons, and Muraki, (1988). This approach is a generalization of the marginal maximum 
likelihood (MML) and Bayes model estimation procedures for unidimensional IRT models 
(see, Bock & Aitkin, 1981, Mislevy, 1986), and has been implemented in TESTFACT (Wilson, 
Wood, and Gibbons, 1991). A comparable model using a logistic rather than a normal-ogive 
^presentation has been studied by Reckase (1985, 1997) and Ackerman (1996a and 1996b). 
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The aim of the present paper is to present a Bayesian approach to estimating normal 
ogive multidimensional IKT models. The method is a generalization of the Bayesian approach 
to estimating the two-parameter normal ogive model (2PNO) by Albert (1992). The method will 
provide marginal posterior density estimates for any parameter of interest. Besides that these 
estimates can be used to judge the accuracy of normal approximations of the standard error of 
MML estimates, the method also has the advantage that restrictions on the dimensionality in the 
MML framework emanating from using Gaussian quadrature do not apply here (Bock, Gibbons, 
& Muraki, 1988). 

This paper is organized as follows. In the next section the procedure will be developed 
for the unidimensional 3-PL normal ogive model (3PNO). In the third section this procedure will 
be generalized to Q-multidimensional normal ogive models, both with and without guessing 
parameters. The former model will be labeled (Q-f 2)PNO, the latter (Q-f- l)PNO. The procedure 
will be defined in such a way that generalization to multiple samples from multiple populations 
is straightforward. In the fifth section, it will be shown how restrictions can be imposed on the 
factor structure and how they can be used to facilitate the interpretation of the model. Then, 
in Section 6, two examples will be given, a simulated example to evaluate parameter recovery 
and an example using real data to exemplify the approach of using restrictions to produce an 
interpretable factor solution. Finally, in the discussion section, some suggestions for further 
research will be given. 



Bayesian Estimation of the 3PNO 

In this paper, a Markov chain Monte Carlo procedure will be used to sample the 
posterior distributions of interest, and the needed chains will be constructed using the Gibbs 
sampler (Gelfand & Smiths, 1990). To implement the Gibbs sampler, the parameter vector 
is divided in a number of components, and each successive component is sampled from its 
conditional distribution given sampled values of all other components. This sampling scheme 
is repeated until the sampled values produce stable posterior distributions. Albert (1992) applies 
Gibbs sampling to estimate the parameters of the well known two-parameter normal ogive 
model (Lawley, 1943, 1944, Lord, 1952, 1953a and 1953b). In this section, the procedure 
will be generalized to the 3PNO, in the following section, this approach will be generalized 
further to a multidimensional framework. In the 3PNO, the probability of a correct response of 
a person i on an item j, denoted by = 1, can be written as 
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P(Yij = = Ij + (1 - 7i)$(»?y) 

= HVij) +7j(l - ®(Vij)), 



( 1 ) 



where 7 ^ is called the ’pseudo-guessing parameter’, 4> denotes the standard normal cumulative 
distribution function, and 77 ^ = a j6{ — /3^ with 6 t the ability of person i , a,, the discrimination 
parameter and (3j the difficulty parameter of the item j, respectively. In (1) the usual expression 
for the 3PNO is rewritten to an expression that supports an interpretation of the model to be 
used below. In this interpretation, there is a probability $>( 77 ^) that the respondent knows the 
item and gives a correct response with probability one, and a probability (1 — <£( 77 ^)) that the 
respondent does not know the item and guesses with 7 ^ as the probability of a correct response. 
So the probability of a correct response is a sum of a term $( 77 ^) and a term 7^(1 - $>( 77 ^)). In 
line with this interpretation, it will prove convenient to introduce a vector of binary variables 
Wij such that 



So if = 0, person i will guess the response to item j> if W{ 3 = 1, person i will know the 
right answer and will give a correct response. Consequently the conditional probability of W \ 3 
given Yij is given by 



To implement the Gibbs sampler the data y, responses of n persons to k items, will be 
augmented with latent data W = ( W n , W nk ) T . Further, following Albert (1992) the data are 
also augmented with latent data Z = (Z 11} ..., Z nfc ) r , where the variables Z i3 are independent 
and normally distributed with mean r\ XJ and a standard deviation equal to one. These variables 
are related to W by Z i3 > 0 if W i3 = 1 and Z i3 < 0 if W i3 = 0. This can be written as 




1 if person i knows the correct answer to item j 
0 if person i doesn’t know the correct answer to item j. 



( 2 ) 



P{Wij = 1 1 Yij = 1 oc 

P(W {j = 0 | = 1 ,^, 7 ,) oc 7; (1 - 
P(W tJ = 1 1 = 0 ,riijCij) = 0 
P(W l] = 0\Y ij = 0,r lij , lj ) =1. 



( 3 ) 



p(Z x] \W tJ>riij ) a 1) (I (Z tj > 0)I(W^ = 1) + I {Zn ^ 0)l(Wij = 0)) (4) 



where Mz, n , v , , ™wmi. estimation MIRT . 5 

evaluated a, Z „ ls assumed ^ ““ % ^ *— - 

«’ «r = MfF, ere given a pnor p® . m. ' " OT " al ' > ' ‘" Stnb ' , " !d ' iKm P™='ers 
parameters are positive, and the pseudo o ' > 0 ' 10 '"“re that the dtscriminatron 

prior Betafl, 1), ^P^esstngparame^^,,,,,^^^^^^ 

The aim of the procedure is to simulate samples from the 
t. 0, Z and W, given by J °' nt P ° Sterior distribution of 



P((,e,z,w M . P(z,w |y; j, 7 , WWJWT) 






( 5 ) 



where p(W^- ) is given by (3) and p(Z I HA _ ) f ,, 

Aithough the distribution g.ven by «, a %) ^ ^ 

distributions of Z, G, £ and 7 are each t t hi ,ntracta ble form, the fully conditional 

conditional d.stribution can be obtatned m five ^ ^ ^ ^ ** fU " 

Step 1 Draw from the conditional distribute of 7 „ 

distribution of Z, cond.t.ona, on W, * and < is given by" ^ ^ ^ ( * ^ 

Z t j |W, 6, £, y dist. I ffiij - 1 ) truncated at the left by 0 ifw . -i 
l truncated at the right by 0 if W t ] = o 

S.ep 2 Draw from ,he condttional disunion of 0 given W 7 „ a „ 

follows that the fully conditional d.s.nbu.ion of e cm , ' ? <4) ““ <5) 11 

- 4 on ay, with , as a tegtesston c “ ^ ^ “ * ** — - 
and cr = l, Therefor, the posKriorore _ . nofma|i that , s _ P"» »* Paramelers p = 0 



< 9 . is distributed n( 

\ 1 At -L 1 2 » 






( 6 ) 



Where 4 ■ (Z « + W2, «? ami variance „i/ £ . 

S-ep 3 Draw from ihe conditional distribution of £ given Rand p Tb 

viewed as coefficients of the regression of Z - = f 7 7 IP ^ He vanables {,. can be 

7 onX=(0- 1))With _ lthe 
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n dimensional column vector with elements -1. Therefore, 

£? I 0 > z j>y is distributed N(£p (X T X) _ 1 )I(a.,- > 0), (7) 

where ^ = (X T X)" 1 X T Z J , 

Step 4 Draw from the distribution of W conditional on all other variables given by (3). 

n 

Step 5 Sample from the conditional distribution of 7 ^ Let tj be defined as tj = £ H W ij = 0) 

i=i 

, that is, tj is the number of persons who do not know the correct answer to item j and guess 
the response. The probability of a correct response of person i on item j given W { j = 0 
is = 1 | Wy = 0) = 7 j . The number of correct responses obtained by guessing, say, 

n 

Sj = Vij> has a binomial distribution, Bin(^, 7 ; ). With the noninformative conjugate 

Beta prior, the posterior distribution of 7 • is 

7 j |W, y - Beta(5 J + 1, tj - Sj + 1). (8) 

The procedure boils down to generating a number of parameter sequences iteratively 
using these five steps. Convergence of the procedure can be evaluated by comparing the 
between and within sequences variance (see, for instance, Gelman, Carlin, Steam & Hall, 
1995). Starting points of the sequences can be provided by the Bayes modal estimates of Bilog 
(Zimowski, Muraki, Mislevy, Sc Bock, 1996). 

Multidimensional normal ogive model 

In a multidimensional normal ogive model with a guessing parameter (Q+2)PNO, the 
probability of a correct response, Yj = 1 , on an item j is given by 



P(Yij = 1 7 j) = 7 j + (1 - 7 



(9) 



o 



where 0* is a vector of Q ability parameters 0 q i, q = 1 = (a 7l , a jq , a jQ , p } ) T 

is the vector of parameters of item j, and r? i; - = a jq@qi ~ Pj- Further, it is assumed that 

all ability parameters i = 1, are independent and have identical multivariate normal 
distributions with mean /j, and covariance matrix £#. 

To identify the model, restrictions have to be imposed on the parameters. One 
approach, which, for instance, is used in the exploratory option of NOHARM (Fraser, 1988), is 
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setting the mean and the covariance matrix equal to zero and the identity matrix, respectively, 

and introducing the constraints a jq = 0, for j = 1 - 1 and q = j + 1 An 

alternative approach is considering the parameters of ability distribution as unknown estimands, 
and to identify the model by setting Q item parameters (3- equal to zero and, for j = 1, ... } Q 

and q = j + 1 , Q , imposing the restrictions a jq = 1 , if j = g, and = 0 , if j ^ q. 

Below it will be shown that these and the restrictions of the exploratory option of NOHARM 
are equivalent. The reason for the latter identification is that it will be easier to generalize the 
estimation procedure to a multiple group situation, where additional ability distributions are 
considered, which are labeled with unique means and covariance matrices. But to keep the 
presentation simple, it will first be assumed that all persons are associated with an identical 
ability distribution; the generalization will be returned to later, together with a generalization to 
incomplete test administration designs, where different persons respond to different items. 

The procedure to be presented will be a generalization of the above procedure. The 
definition of the variables, W = {W lu W nk ) and Z = ( Z n , ..., Z nk ) is as above, but with 
an alternative definition of 77 ^. It is assumed that the ability of person n, 0 n , is Q-variate 
normally distributed with mean (i and variance £, that is, 6 n ^ N(fj, , S 0 ). Further, the prior 
for the item parameters £ is given by 7r(£) = IlljLi I {<*jq > 0), that is, it is assumed 

that the probability of a correct response is positively related to all latent abilities involved. 
As above, the pseudo guessing parameter 7 ^ will have a Beta prior. The joint posterior of 
£, 0 ,Z,W, 7 is given by 



Pit 0 , Z, W, 7 , /*, E |y ) = p( Z, W |y ; 7, *,)p(* I M, |S, )p(E»)p(0?(7) 

= n n p( z h \W n ; Vij)p(w tj \ yij \ r] iv 'y j )p{e t |/x, s e ) 

t— 1 JS=1 

p(m|e* )p(S®)p(€)p(7) 

where p{Z {j \Wij ; 77 ^) and p{Wij \yij ; 77 ^, 7 ^) are as defined in (4) and ( 5 ), with the proper 
substitution of 77 ^. Finally, p(pt |S^ ) and p(S^) will be normal and inverse Wishart, 
respectively. 

Again, the procedure consists of simulating from fully conditional posteriors, in this 
case, the fully posterior distributions of £, 6 , Z, W, 7 , //, Therefore, the procedure consists 

of seven steps: 

(1) Draw conditional on 6. 

(2) Draw \x conditional on £0 and 6 

ERIC 
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(3) Draw W conditional on 6 , £, 7 and Y. 

(4) Draw Z conditional on W, 6 , and £. 

(5) Draw 0 conditional on Z, £, E#, and pi. 

( 6 ) Draw f conditional on Z and 

(7) Draw 7 conditional on W and Y. 

In this section, it will be assumed that there are appropriate initial estimates of 0, £ and 
7 , how these initial estimates might be obtained will be discussed later With this assumption, 
the seven steps can be executed in the following order 

Step 1 The conjugate prior distribution for (/x, E 0 ) is a product of a normal and an inverse- 
Wishart distribution (see, for instance, Box & Tiao, 1973). Let A 0 be the scale matrix and let 
vq be the degrees of freedom for the prior of E^. Further, let /x 0 be the prior mean and let ft 0 be 
the number of prior measurements on the E$ scale. The posterior is now also inverse-Wishart 
distributed with parameters 



Mn 


— *o n .. _l_ 

Ko+n^O 


Vn 


= vo + n 


Kn 


= fto + to 


A n 


= A 0 + S 4- 



*o+n 



Ko+n 



(0 - n 0 ){0 - v 0 ) T , 



n __ __ 

where S is a matrix with elements = ^2{9 ip - 9 p )(9 iq — 9 q ) and 9 P and 9 q are the means 

1=1 

of the ability parameters of dimension p and q, respectively. The choice of a value for fi 0 
will be returned to below. With respect to the choice of k 0 ,v 0 and A 0 , a noninformative prior 
distribution is obtained if ft 0 — ► 0,vo — ► “1 and |A 0 | — ► 0. This results in the multivariate 
version of Jeffreys’ prior density. The corresponding posterior distribution can be written as 



E 0 10 ~ Inverse-Wishart^! (S) 

~ N(d,J2 g /n). 



( 10 ) 



So using Jeffreys’ prior E^ can be drawn from the Inverse-Wishart distribution. 

Step 2 Draw pi from its full conditional distribution. Using Jeffreys’ prior results in /x | E^, y ~ 

N(d,H e /n). 



Step 3 Sampling W from its full conditional distribution proceeds as in Step 4 of the previous 
section, with the proper substitution of 77 ^. 
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Step 4 Sampling Z from its full conditional distribution proceeds as in Step 1 of the previous 
section, with the proper substitution of 

Step 5 To draw from the conditional distribution of 0 , an orthogonal standardized ability 0° is 
defined. So the elements of 0° = ( 6 ° 1 } ... 0 °q) t have independent standard normal distributions. 
Let L be the Cholesky decomposition of £#, that is, £$ = LL T . Define the orthogonal ability 
vector for person i , 6° = pi). Now can be written as 

Q 

Vij = 'y 'j SXjqQgj) ~ Pj 

9=1 

Q Q 

= LhqO q j + /ig) — Pj 

q=l h=l 

and in matrix notation as 



r) x = ALlT^-pi + /i) - / 3 = A(L0? + pi) - /3, 

with 77^ and /3 vectors of length 0*, 0° and pi vectors of length Q and A a k x Q matrix with 
elements Otj q . The ability parameters 6° has a full posterior density given by 

k 

p(0i |£, Z, y, W) cx 4 >( 0 °\ o, I) n HZijWij, 1). 

i= i 

This can be written as 



Zj + /3 — Apt — B0° + £*, 

where B = AL and e* is a vectors elements which are iid iV(0, l).It then follows that 

0° is distributed N ((I + £- 1 )- 1 ir 1 0°, (I + IT 1 )" 1 ) , 

with 0° the common least squares estimate 0° = (B T B) _1 B T (Z f + /3— Apt) and £ = 
(B r B) -1 . Now 6 { can be obtained by the transformation L0^ + pi = 0 t . 

Step 6 Conditional on Z and 0, the posterior distributions of are independent with 

EKLC 
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n Q Q 

p(€j \W y 0 , Z, y) oc Y[ ; ^2 a 3 q^iq ~ l) fj I(Qj g > 0). 

i=rl q—l q ~ 1 

Define X = (0 1} 0 q , O q , —1) with = (0i 9 , ...# n9 ) r and -1 is an n dimensional column 
vector with elements -1. Conditional on 9 , Z, = (Z y , Z nj -) T satisfies the linear model 

Zj = + Ejy 

where £j = ( £ljl ...,e flJ -) T and the e*j are lid AT(0, 1). It then follows that 

Q 

\°r Z > y is distributed iV(^, (X^X)" 1 ) JJ I(a,- 9 > 0), 

9=1 

where = (X T X) -1 X T Zj. 

Step 7 Sampling 7^ from its full conditional distribution proceeds as in Step 5 of the previous 
section. 



As above, these steps can be used to iteratively generate a number of parameter 
sequences and convergence can be evaluated by comparing the between and within sequences 
variance. Finding starting values for the procedure may be obtained from NOHARM or 
TESTFACT, but an alternative approach will be sketched in the next section. 

Above, it was already mentioned that there are several approaches to identify the 
model. Consider the model with ability parameters 9 t and identifying restrictions fj, = 0 , 
a 3Q = !. 3 = 9, and a jq = 0, if j ± q , for j = 1, ..., Q and q = j + 1 , ..., Q, on one 

hand, and the model with orthogonal ability parameters 9° and identifying restrictions jj. = 0 , 
£0 — I, and = 0, for j = 1, ..., Q — l and <7 = j ' + 1 , Q t on the other hand. Let A and 
A be the matrices of factor loadings for the former and latter parameterization, respectively. 
Analogously to step 5, 9 Z can be transformed to 9° by 9° = L ~ l 9 u where L is the Cholesky 
decomposition of £#. Because Lis lower triangular and A9 Z = AL 9° = A°9°, the restrictions 
a jq ~ 1) if j — 9> an d otj q = 0, if J i= q, for j = 1 , ..., Q and q = j + 1, ..., Q, are transformed 
into restrictions a° jq = 0, for j = 1 , .... Q - 1 and q = j + 1 , ..., Q . On the other hand, 
defining the lower triangular matrix F as the first Q rows of A° and applying 0 * = F9°, results 
in £ 0 = F T F 1 and A = A°F _1 , which in turn produces restrictions a j q = 1 , if j = q , and 

1 9 
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ajq = °’ ' {j * q ' f ° r j = h •’ Q and * = 3 + 1, .... Q. So the two parameterization of the 
model are easily interchanged. 

As was already ment.oned above, the procedure easily generalizes to the case of 
multiple groups and incomplete designs. In the case of multiple groups, the mean and the 
covariance matrix of the ability distribution of one of the groups can be set equal to zero 
and the identity matrix, respectively. So for this group the first two steps are not performed, 
these steps are only performed to generate means and covariance matrices for the remaining 
groups. Restricting the mean and the covariance of one of the ability distributions suffices for 
f.xmg the location and the scale of the latent continuum. However, this is not sufficient to 
identify the model, because the model probabilities (1) are left unaltered under a permutation 
of the ability dimensions. Therefore, additional constraints a jq = O.for j = l and 

q = j + l , .... Q are needed to completely identify the model. Steps 3 and 4 are performed per 
item, so generalization to multiple groups and incomplete designs only entails that the variables 
Wij and Zij P erformed *e cases where person * responded to item j. Step 5 is performed 
per person, so here both the group to which the person belongs and the items made should be 
taken into account. Finally, the steps 6 and 7 for generating item parameters should take into 
account the persons actually confronted with the items. 



Subscale Factor Analysis 



In the previous section it was aigued that the model can be identified by setting a mean 
and a covariance matnx equal to zero and the identity matrix, respectively, and introducing the 
constraints a jq - 0,for j = 1, 1 and? = j+ 1 , These are the constraints also used 

in the exploratory option of NOHARM (Fraser, 1988). In this approach, the mode, is identified 
by assuming that the responses on the first item are uniquely determined by the first ability 
dimension, the responses on the second items are uniquely determined by a mixture of the first 
and the second ability dimension, and so forth. In general, these identification restriction will be 
of little help in providing an interrelation for the ability dimensions. Therefore, in a exploratory 
factor analysis, the factor solution is usually visually or analytically rotated. Often, the rotation 
scheme is devised to approximate Thurstone’s simple structure criterium (Thurstone, 1947), 
where the factor loadings are split into two groups, the elements of the one tending to zero 
and the elements of the other toward unity. The approach suggested in this section has much 
in common with this approach: the idea is to identify the dimensions with subscales of items 
loading on one dimension only. The idea is to either identify these S <Q subscales on a pnon 
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grounds, or to identify these subscales using an iterative search based on fitting a unidimensional 
IRT model by discarding items from the scale. In the latter case, there are two reasons for 
choosing the Rasch model as the unidimensional IRT model. Firstly, the theory evaluation of 
fit to the Rasch model is well founded and the asymptotic distributions of the many fit statistics 
are well known (Andersen, 1973, Martin Lof, 1973, 1974, Kelderman, 1984, 1989, Molenaar, 
1983, Glas, 1988, Glas & \fcrhelst, 1989, 1995). However, this argument is not completely 
overwhelming, because procedures for evaluation of model fit are also available in the general 
MML framework (Yen, 1981, 1984, Mislevy & Bock, 1990, Reiser, 1996, Glas, 1998) and 
the framework of non-parametric IRT (Rosenbaum, 1984, Holland & Rosenbaum, 1986, Stout, 
1987, 1990 and Junker, 1991). A second and more important reason relates to the intuitive 
concept of a measurement scale: a continuum where the order of the outcomes of a measurement 
does not depend on the objects positions on the scale. Translated to a test with dichtomous items: 
the probabilities of a correct response on items referring to a scale is independent of the position 
on the scale. This intuitive requirement is formalized by Rasch (1967) and leads to the Rasch 
model (Fischer, 1995). 

These considerations lead to the following heuristic: 

(1) Set counter 5 = 1 

(2) Estimate the parameters of the items in the target set; 

(3) Compute fit statistics and remove poorly fitting items from the target set to the set ” removed 
items”; 

(4) If the overall model fit is unacceptable, goto 2; 

(5) If 5 = 5 then stop else ” removed items” becomes the target set, s = s + 1 and goto 2. 

In the first iteration of the selection process ( s = 1), the target distribution can 
be set to the complete test, and the inner loop (steps (2)-(4)) is repeated until an acceptable 
unidimensional subscale is found. In the next iteration (5 = 2), the nest subscale is selected 
from the remaining items, and so forth, until 5 = 5. Detailing this procedure further, the item 
parameters estimates of step (2) can, for instance, be computed using conditional maximum 
likelihood (see, for instance, Molenaar, 1995), which has the advantage that these estimates 
do not depend on an assumption with respect to the ability distribution. These estimates can 
be computed using RSP (Glas & Ellis, 1993), which program also computes the so-called R\ 
and Rq statistics (Glas, 1988) that can be used in step (3) and (4). The i?i-test can be used for 
evaluation of the form of the ICC’s and the J? 2 -test for evaluation of local independence. 



After identification of 5 sets of scaled items, there will usually be a set of remaining 
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items that may load an all ability dimensions. Of course, the number of items in this set depends 
on the size of the critical regions used in step (3) and (4). Next, restrictions will be imposed on 
the matrix of factor loadings A to reflect the structure of the subscales found, that is, if item j 
belongs to subscale q, ctj q = 1 , and a.j q , — 0, for q f = 1 , Q } q f q. Finally the remaining 
dimensions q = S + 1 , , ..., Q are identified using ay q = l,and a ; y = 0 , for some item f not 
belonging to a subscale, and q f = 1 , ..., Q , q f ^ q. 

Examples 

In this section, two studies of the performance of the MCMC method will be presented, 
one for the case of the 3PNO and one for the case of the (3+2)PNO. For the first study, data 
were simulated using the 3PNO parameter values given in Table 1 . These parameters are a 
subset of the items parameters of a real, and therefore realistic data set (ACT, 1997). Using 
these parameters, ability parameters for n = 1000 simulees were drawn from a standard normal 
ability distribution and for every simulee a response patterns was generated using the 3PNO. 
Next, the items parameters were estimated using the MCMC procedure described above and 
Bilog-MG. In both procedures the same Beta(a, b) prior on 7 was used. The values for a and 
b were set equal to 20 * y true + 1 and 20 * (1 —' 7 true) + respectively, the values 7 true are 
listed in Table 1. Further, in the Bilog analysis a log normal (0,0.5) prior on a was used. The 
MCMC procedure had a run length of 30000 iterations. The starting values used in the MCMC 
procedure were, a* = 1, 0 { = 0, 7 * = 7 itrue for all items i, and 6j = Ofor all persons j. In areal 
estimation situation those starting values can be improved by the use of, for example, Bilog- 
MG estimated parameter values. From examining the plots of sampled parameter values, it was 
concluded that a bum-in period of 1000 iterations was sufficient. These first 1000 iterations 
were discarded. 



Insert Table 1 about here 



The mean of the generated posterior distribution and posterior standard deviation of 
the parameters issued from the MCMC procedure are given in Table 2. In Table 3, the MML 
estimated parameter values and their standard errors issued from Bilog-MG are given. The 0 
parameters in Table 3 are transformed 0 = a/? bjlog , because in Bilog-MG the parameterization 
$(a(0 — 0)) is used. Comparing the Tables 1, 2 and 3, it can be seen that parameter recovery 
of the MCMC and Bilog-MG procedures are comparable. 



MCMC estimation MTRT 



14 



.. eva Iuate (he performance of the MCMC 

imensional ( 3 +2 )PNO model, 20 0 0 response PrOC ' dUr ' ,n “^ma.ion „, th the 3. 
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° As wilh "» 3PNO model, from examina, T “* P ' rf ° m ' ed 1 ™ b a « lengih 
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In Table 5 the results of the MPlur 
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msnlts in Table , ,„ d Tab, e 7 ,Prom mspecdon ^ “ "* *" •» 0. 
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Table 1 



True item parameter values 



ITEM 


a 


p 


7 


1 


0.642 


-1.619 


0.19 


2 


0.806 


-1.533 


0.15 


3 


0.956 


-1.292 


0.11 


4 


0.972 


-1.061 


0.14 


5 


1.045 


-0.245 


0.37 


6 


0.834 


-0.264 


0.14 


7 


0.614 


0.023 


0.17 


8 


0.796 


0.213 


0.10 


9 


1.171 


-0.669 


0.19 


10 


1.514 


0.480 


0.31 
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Table 2 

Recovery of 3PNO parameters by MCMC 



ITEM 


a 


P 


7 


er ( a ) 


<*( 0 ) 


< rb ) ■ 


1 


0.564 


- 1.607 


0.220 


0.111 


0.114 


0.087 


2 


0.818 


- 1.500 


0.187 


0.125 


0.119 


0.081 


3 


0.928 


- 1.185 


0.161 


0.142 


0.112 


0.079 


4 


0.826 


- 0.829 


0.173 


0.125 


0.111 


0.080 


5 


0.927 


- 0.399 


0.325 


0.202 


0.166 


0.079 


6 


0.975 


- 0.369 


0.113 


0.130 


0.095 


0.049 


7 


0.629 


0.038 


0.214 


0.151 


0.200 


0.082 


8 


1.374 


0.682 


0.221 


0.444 


0.332 


0.048 


9 


1.093 


- 0.691 


0.160 


0.173 


0.107 


0.061 


10 


1.027 


0.273 


0.279 


0.234 


0.198 


0.057 
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Recovery of 3PNO 



ITEM a 



1 

2 

3 

4 

5 

6 

7 

8 

9 

10 



0 



Table 3 

parameters by BILOG 



MG 



0.653 
0.865 
0.919 
0.809 
, 0.879 
0.930 
0.582 
1.139 
1.029 
0.942 



-1.736 

-1.600 

-1.249 

-0.869 

-0.445 

-0.416 

-0.024 

0.499 

-0.730 

0.203 



7 c{a) a(p) g(y) 



0.191 

0.150 

0.119 

0.145 

0.311 

0.084 

0.195 

0.199 

0.131 

0.267 



0.113 

0.125 

0.122 

0.108 

0.152 

0.116 

0.106 

0.279 

0.135 

0.209 



0.165 

0.161 

0.143 

0.121 

0.172 

0.099 

0.088 

0.150 

0.140 

0.150 



0.088 

0.080 

0.074 

0.077 

0.082 

0.047 

0.081 

0.049 

0.061 

0.063 
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ITEM 



Table 4 

True item parameters Q+21PNn 



a i 



a 2 






0 



2 

0.19 
0.15 
0.11 
0.14 
0.37 
0.14 
0.17 
0.10 
0.19 
0.31 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 
0.20 



1 

2 


1.000 

.000 


.000 

1.000 


3 


.000 


.000 


4 


1.000 


.000 


5 


1.000 


.000 


6 


1.000 


.000 


7 


1.000 


.000 


8 


1.000 


.000 


9 


1.000 


.000 


10 


1.000 


.000 


11 


1.000 


.000 


12 


.000 


1.000 


13 


.000 


1.000 


14 


.000 


1.000 


15 


.000 


1.000 


16 


.000 


1.000 


17 


.000 


1.000 


18 


.000 


1.000 


19 


.000 


1.000 


20 


.000 


.000 


21 


.000 


.000 ; 


22 


.000 


.000 J 


23 


.000 


.000 ] 


24 


.000 


.000 1 


25 


.000 


.000 1 


26 


.200 


.800 


27 


.200 


.200 


28 


.800 


.400 



1 .000 
.000 
1.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
.000 
1.000 
1.000 
1.000 



.000 

.600 

.200 



.000 
.000 
.000 
.500 
1.000 
- 1.000 
-1.000 
.500 
1.000 
-1.000 
-1.000 
.500 
.250 
-.250 
-.500 
.500 
.250 
-.250 
-.500 
1.000 
.000 
1.000 
1.000 
.000 
1.000 
.000 
1.000 
1.000 
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ITEM 

1 

2 

3 

4 

5 

6 

7 

8 

9 

10 
11 
■12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 

27 

28 



Table 5 

Recovery of 3PNO parameters by MCMC procedure with prior on a 



Ot 1 


a 2 


a 3 


p 


7 


1.148 


0.000 


0.000 


0.023 


0.188 


0.138 


0.927 


0.000 


-0.011 


0.150 


-0.027 


0.046 


1.007 


0.046 


0.111 


0.919 


0.034 


0.046 


0.504 


0.144 


0.854 


0.088 


0.110 


0.923 


0.378 


0.881 


0.048 


0.029 


-0.923 


0.142 


1.061 


0.040 


0.062 


-1.072 


0.172 


0.989 


0.066 


0.028 


0.463 


0.105 


1.272 


0.075 


0.129 


1.258 


0.190 


0.982 


0.031 


0.057 


-1.067 


0.313 


0.977 


0.047 


0.102 


-1.139 


0.201 


0.034 


1.083 


0.132 


0.652 


0.202 


0.040 


1.015 


0.102 


0.294 


0.202 


0.043 


1.034 


0.035 


-0.270 


0.203 


0.078 


0.943 


0.050 


-0.477 


0.201 


0.094 


1.088 


0.130 


0.556 


0.203 


0.038 


1.089 


0.180 


0.334 


0.205 


0.061 


0.984 


0.076 


-0.226 


0.203 


0.035 


1.020 


0.028 


-0.459 


0.202 


0.062 


0.102 


1.025 


L041 


0.206 


0.078 


0.054 


1.330 


0.039 


0.203 


0.039 


0.062 


1.135 


1.146 


0.212 


0.056 


0.086 


1.207 


1.221 


0.209 


0.105 


0.040 


1.032 


0.009 


0.204 


0.070 


0.061 


1.356 


1.355 


0.213 


0.253 


0.782 


0.040 


-0.015 


0.202 


0.096 


0.182 


0.607 


1.069 


0.200 


0.716 


0.391 


0.201 


-0.933 


0.201 



<7(0:1) 


o-(a 2 ) 


o-(a 3 ) 


a (P) 


<7(7) 


0.109 






0.059 


0.012 


0.050 


0.074 




0.049 


0.011 


0.053 


0.051 


0.090 


0.050. 


0.010 


0.093 


0.027 


0.035 


0.068 


0.011 


0.162 


0.061 


0.072 


0.165 


0.015 


0.079 


0.034 


0.025 


0.055 


0.011 


0.098 


0.032 


0.043 


0.067 


0.012 


0.098 


0.042 


0.024 


0.066 


0.010 


0.213 


0.053 


0.074 


0.197 


0.011 


0.110 


0.027 


0.042 


0.073 


0.015 


0.093 


0.036 


0.055 


0.068 


0.013 


0.029 


0.109 


0.064 


0.086 


0.012 


0.032 


0.096 


0.055 


0.066 


0.013 


0.033 


0.088 


0.029 


0.051 


0.013 


0.045 


0.080 


0.037 


0.051 


0.013 


0.053 


0.114 


0.062 


0.086 


0.012 


6.031 


0.106 


0.065 


0.073 


0.013 


0.041 


0.081 


0.047 


0.052 


0.013 


0.028 


0.085 


0.024 


0.051 


0.013 


0.047 


0.062 


0.141 


0.139 


0.013 


0.051 


0.041 


0.143 


0.068 


0.013 


0.035 


0.048 


0.209 


0.208 


0.014 


0.046 


0.059 


0.207 


0.203 ■ 


0.013 


0.053 


0.032 


0.102 


0.058 


0.013 


0.055 


0.049 


0.286. 


0.276 


0.013 


0.054 


0.071 


0.031 


0:051 


0.013 


0.057 


0.068 


0.089 


0.102 


0.012 


0.069 


0.057 


0.056 


0.053 


0.013 
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4 

5 

6 

7 

8 

9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 
21 
22 

23 

24 

25 

26 
27 

28 



gi 

1.000 

0.000 
0.000 

0.923 
0.890 
0.896 
1.071 
0.991 
1.310 
1.020 
0.985 
- 0.182 
- 0/130 
- 0.124 
- 0.044 
- 0.047 
- 0.156 
- 0.076 
- 0.145 
0.024 
0.082 
- 0.052 
- 0.004 
0.125 
0.016 
0.159 
0.066 
0.692 

Covariance 
1 X 975 
0.016 0.938 



02 

0.000 
1.000 
0.000 
0.025 
0.104 
0.064 
0.043 
0.103 
0.104 
- 0.010 
0.052 
1.137 
1.036 
1.059 
0.967 
1.113 
1.111 
1.018 
1.058 
0.052 
- 0.018 
- 0.041 
0.007 
- 0.037 
- 0.062 
0.822 
0.166 
0.433 



0^3 

0.000 
0.000 
1.000 
- 0.042 
0.056 
- 0.093 
- 0.017 
- 0.098 
0.066 
- 0.035 
0.046 
0.092 
0.070 
- 0.059 
- 0.009 
0.096 
0.152 
0.034 
- 0.087 
1.030 
1.341 
1.139 
1.226 
1.038 
1.320 
- 0.033 
0.599 
0.159 



0 

0 . 004 ' 
- 0.006 
0.045 
0.506 
0.956 
- 0.930 
- 1.071 
0.465 
1.287 
- 1.079 
- 1.138 
0.673 
0.294 
- 0.271 
- 0.474 
0.552 
0.333 
- 0.224 
- 0.458 
1.042 
0.037 
1.132 
1.228 
0.010 
1.310 
- 0.014 
1.062 
- 0.938 



7 

0.184 
0.152 
0 . 11 ] 
0.144 
0.380 
0.142 
0.172 
0.105 
0.190 
0.313 
0.201 
0.204 
0.202 
0.202 
0.202 
0.203 
0.205 
0.203 
0.203 
0.206 
0.203 
0.211 
0.209 
0.204 
0.212 
0.202 
0.200 
0.201 



0.100 0.064 0.067 

0 186 0.095 0 098 

0084 0.062 0.064 

0I °5 0.070 0.073 
0105 0.065 0.067 

0-238 0.096 0 099 

0119 0.072 0.077 

0.098 0.069 0.072 
0.082 0.134 0.080 

0069 0.100 0.071 

n!S 0095 a07 ° 

nnS 0 084 0 066 

0 074 0.120 0.078 

OnS 0 119 0078 

0 064 0.092 0.069 

0069 0.095 0.071 

0080 0.084 0 151 

0079 0.081 0.154 

onS 0086 0 166 

0.094 0.099 0.204 

0.067 0.069 0.109 

0098 0.102 0 193 

0061 0.077 0.064 

0 071 0.074 0.091 

0 £ 75 __ 0066 _ 0064 ( 



0.052 
0.049 
0.049 
0.068 
0.182 
0.054 
0.068 
0.064 
0.210 
0.074 
0.068 
0.099 
0.066 
0.052 
0.050 
0.084 
0.073 
0.053 
0.052 
0.144 
0.069 
0.159 
0.198 
0.059 
0.183 
0.051 
0.104 
0.053 



0.012 
0.011 
0.010 
0.011 
0.015 
0.011 
0.012 
0.010 
0.012 
0.015 
0.012 
0.012 
0.012 
0.013 
0.013 
0.012 
0.013 
0.013 
0.013 
0.013 
0.013 
0.012 
0.013 
0.013 
0.012 
0.013 
0.012 
0.013 
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Table 7 

NOHARM estimates (3+2)PNO 



ITEM 


Oil 




<*3 


0 


1 


1.000 • 


0.000 


0.000 


-0.038 


2 


0.000 


1.000 


0.000 


0.014 


3 


0.000 


0.000 


1.000 


-0.042 


4 


0.686 


0.007 


-0.014 


-0.480 


5 


0.577 


0.075 


-0.006 


-0.801 


6 


0.671 


0.042 


-0.056 


0.922 


7 


0.783 


0.015 


0.007 


1.049 


8 


0.734 


0.064 


-0.081 


-0.436 


9 


0.958 


0.021 


0.130 


-1.230 


10 


0.686 


-0.040 


0.036 


1.025 


11 


0.765 


0.010 


0.081 


1.143 


12 


-0.117 


1.128 


0.113 


-0.634 


13 


-0.073 


1.068 


0.023 


-0.284 


14 


-0.114 


1.070 


-0.061 


0.280 


15 


-0.042 


0.978 


-0.032 


0.478 


16 


-0.037 


1.108 


0.066 


-0.526 


17 


-0.096 


1.082 


0.125 


-0.297 


18 


-0.069 


1.033 


0.032 


0.237 


19 


-0.107 


1.055 


-0.085 


0.463 


20 


0.032 


0.019 


0.956 


-0.962 


21 


0.072 


-0.070 


1.293 


-0.020 


22 


-0.040 


-0.056 


0.922 


-0.942 


23 


-0.033 


0.020 


1.061 


-1.066 


24 


0.098 


-0.062 


0.968 


0.007 


25 


0.014 


-0.024 


1.022 


-1.043 


26 


0.127 


0.818 


-0.019 


0.024 


27 


0.065 


0.161 


0.645 


-1.081 


28 


0.498 


0.389 


0.195 


0.918 


Covariance 


1.638 










0.062 


0.886 








-0.070 


0.042 


0.995 
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Table 8 



TESTFACT estimates (3+2)PNO 



ITEM 


Oi 1 


a 2 




p 


1 


1.000 


0.000 


0.000 


0.013 


2 


0.000 


1.000 


0.000 


-0.013 


3 


0.000 


0.000 


1.000 


0.046 


4 


0.782 


-0.004 


-0.041 


0.467 


5 


0.678 


0.082 


0.022 


0.809 


6 


0.760 


0.029 


-0.083 


-0.928 


7 


0.883 


-0.005 


-0.021 


-1.051 


8 


0.840 


0.071 


-0.106 


0.429 


9 


1.035 


0.033 


0.073 


1.158 


10 


0.880 


-0.060 


-0.027 


-1.086 


11 


0.818 


-0.005 


0.042 


-1.116 


12 


-0.181 


1.203 


0.097 


0.651 


13 


-0.131 


1.084 


0.040 


0.283 


14 


-0.132 


1.099 


-0.069 


-0.271 


15 


-0.066 


1.014 


-0.010 


-0.470 


16 


-0.061 


1.146 


0.112 


0.529 


17 


-0.156 


1.130 


0.142 


0.307 


18 


-0.085 


1.024 


0.023 


-0.225 


19 


-0.158 


1.097 


-0.090 


-0.455 


20 


0.037 


0.063 


1.008 


0.991 


21 


0.095 


-0.032 


1.251 


0.022 


22 


-0.033 


-0.069 


1.047 


1.011 


23 


0.004 


0.016 


1.031 


1.035 


24 


0.135 


-0.053 


0.972 


-0.004 


25 


0.046 


-0.061 


1.060 


1.066 


26 


0.109 


0.836 


-0.041 


-0.025 


27 


0.076 


0.142 


0.587 


1.039 


28 


0.572 


0.409 


0.155 


-0.924 


Covariance 


1.175 










0.081 


0.754 








-0.032 


0.037 


0.890 
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